Characterizing Spinning Black Hole Binaries in Eccentric Orbits with LISA 
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The Laser Interferometer Space Antenna (LISA) is designed to detect gravitational wave signals 
from astrophysical sources, including those from coalescing binary systems of compact objects such 
as black holes. Colliding galaxies have central black holes that sink to the center of the merged galaxy 
and begin to orbit one another and emit gravitational waves. Some galaxy evolution models predict 
that the binary black hole system will enter the LISA band with significant orbital eccentricity, 
while other models suggest that the orbits will already have circularized. Using a full seventeen 
parameter waveform model that includes the efi'ects of orbital eccentricity, spin precession and 
higher harmonics, we investigate how well the source parameters can be inferred from simulated 
LISA data. Defining the reference eccentricity as the value one year before merger, we find that 
for typical LISA sources, it will be possible to measure the eccentricity to an accuracy of parts in 
a thousand. The accuracy with which the eccentricity can be measured depends only very weakly 
on the eccentricity, making it possible to distinguish circular orbits from those with very small 
eccentricities. LISA measurements of the orbital eccentricity can help constraints theories of galaxy 
mergers in the early universe. Failing to account for the eccentricity in the waveform modeling can 
lead to a loss of signal power and bias the estimation of parameters such as the black hole masses 
and spins. 



I. INTRODUCTION 



Binary systems of compact objects will be ubiqui- 
tous sources for the Laser Interferometer Space Antenna 
(LISA) P, Q . Observations have shown that today there 
are massive black holes in the center of nearly all galax- 
ies 0-lH . When galaxies collide their central black holes 
sink to the center of the merged galaxy and begin to orbit 
one another, losing energy and angular momentum in the 
form of gravitational waves 0] . Gravitational wave (GW) 
emission rapidly erases any initial eccentricity 0, Q i so it 
has long been thought that eccentricity could be ignored 
when modeling the signals from massive black hole bi- 
naries 0. More recently, however, it has been shown 
that the mechanisms that may harden the binary to the 
point where the gravitational wave emission takes over 
all tend to drive up the eccentricity [l0l - [l9| . The ques- 
tion then is whether significant eccentricity survives until 
the final year or so before merger. Figure [1] shows the 
eccentricity evolution as a function of orbital frequency 
for systems that enter the gravitational wave dominated 
evolution stage at frequency /gw with eccentricities of 
Bgw = 0.95 and Cgw = 0.5. The tracks are computed us- 
ing the leading-order Peters and Matthews 0] evolution 
equations. These equations predict that once the eccen- 
tricity drops below e ~ 0.3, it decays as e ~ j-i9/i8 _ 
other words, roughly a factor of ten in eccentricity is lost 
per decade of frequency. The rate of decay is slower for 
systems with very high eccentricities, allowing them to 
maintain significant eccentricity for longer. For typical 
LISA black hole binaries, the GW decay time drops be- 
low the Hubble time when the systems are 3 to 5 decades 
in frequency from the LISA band, so unless egw is very 
close to unity, a purely GW driven orbital evolution re- 
sults in nearly circular orbits in the LISA band. On the 
other hand, if the hardening mechanism (e.g. gas dy- 



namics or stellar scattering) continues to dominate the 
dynamics until the system is close to the LISA band, 
then interesting eccentricities can be maintained. In a 
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FIG. 1: Gravitational wave driven eccentricity evolution as a 
function of orbital frequency. 

recent study [l^, Sesana has shown that stellar scat- 
tering produces LISA sources with eccentricities in the 
range eo = 10~^ — ?> 0.2. Moreover, as shown in Figure 
8 of Ref. [l^l , the distribution of eccentricities depends 
on the component masses in a particular way, so measur- 
ing this distribution can help constrain black hole merger 
models. While the distribution of component masses and 
spins as a function of luminosity distance will likely play 
a more important role in studying galaxy - black hole 
co-evolution j2l| - |23j . the eccentricity distribution may 
provide useful additional constraints. These considera- 
tions suggest that it is desirable to include the effects 
of eccentricity in the gravitational waveform, bringing 
the total number of parameters needed to describe a 
black hole inspiral to seventeen: the redshifted black hole 
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masses (mi, 7712), the distance to the source (-Dl), the ini- 
tial radial eccentricity and semi-major axis (eo,ao), the 
dimensionless spin parameters for the two black holes 
(XijXa); the source sky location (cos0, </>), the initial 
orientation of the angular momentum and spin vectors 
(cos , ) , (cos 9si , (f>Si ) , (cos 7 </'S2 ) I sjnd two initial 
phase parameters (no,0o)- Failing to include the eccen- 
tricity in the waveform model will lead to a loss of signal- 
to- noise [U, 25 1 and to biases in the recovery of the other 
parameters 26[ [27|. 

The instantaneous gravitational waveforms describing 
the inspiral of a black hole binary with arbitrary spins, 
masses, and orbital eccentricity were calculated by Kid- 
der [23 to second post-Newtonian order (2 PN, or order 
v'^/c'^ in the relative velocity), and extended to 2.5 PN 
order by Faye, Blanchet and Buonanno [l^H^]- Majar 
and Vasiith [31[ later introduced a convenient frame- 
work for expressing the waveforms in a form amenable 
to producing detection templates. Together with a so- 
lution to the conservative equations of motion, adiabati- 
cally advanced with the angular momentum and energy 
dissipation equations [l^, we can build the time depen- 
dent gravitational waveforms for a general binary black 
hole system with the full seventeen parameters neces- 
sary to describe the system. LISA observations of binary 
black hole inspirals can be compared to these waveforms 
to produce posterior distributions for the model parame- 
ters. These observations will allow us to constrain galaxy 
merger scenarios [2l| - [23| . and allow us to test General 
Relativity in the dynamical strong field regime [s^l ■ Left 
unaccounted for, the effects of orbital eccentricity may be 
mistaken for a departure from General Relativity (in par- 
ticular, even small eccentricities that produce negligible 
power in higher harmonics can lead to easily detectable 
changes in the phase evolution of the signal). 

Black hole binary systems in eccentric orbits may also 
be detected by ground based gravitational wave detec- 
tors such as the Laser Interferometer Gravitational wave 
Observatory (LIGO) and Virgo [s^. Some models even 
predict that inspiral signals may enter the LIGO band 
with e > 0.9 and that eccentric templates will be nec- 
essary to detect such sources Current LIGO data 
analysis uses circular templates and may need to be gen- 
eralized to include eccentricity. 

We have previously described a method for combin- 
ing the instantaneous gravitational waveforms for eccen- 
tric binary systems with a post-Keplerian solution to the 
equations of motion that is adiabatically advanced us- 
ing the the orbit-averaged dissipation and spin precession 
equations to build ready-to-use gravitational waveforms 
for the general case of a spinning black hole binary system 
in an eccentric orbit [35j . Here we present the results of a 
parameter estimation study for spinning eccentric binary 
black hole sources for the LISA mission. This is an exten- 
sion of the Lang and Hughes LISA parameter estimation 
study of spinning binary black holes in circular orbits [36j 
and is the first to include the full seventeen parameters 
that describe a general black hole binary inspiral. 



II. WAVEFORM MODEL 

The equations needed to numerically calculate time de- 
pendent gravitational waveforms for binary black hole 
systems have been computed to 2.5PN order in the am- 
plitude and phase [28l - l3]| . The resulting system of equa- 
tions could be numerically evolved to produce waveforms 
for our study, but the computational cost of resolving the 
motion of the black holes on the orbital timescale is too 
large for the parameter estimation study we wish to per- 
form. Taking advantage of the separation of time scales 
in the waveform model we have developed an efficient 
waveform generator at 1.5 PN order in the amplitude 
and phase [35|- These waveforms include the effects of 
periastron precession, the precession of the orbital plane 
due to spin-orbit coupling, and higher harmonics from 
the higher order mass and current multipole moments. 
Our waveform model does not include the effects of spin- 
spin coupling which enter at 2PN order. In future work 
we plan to extend the waveform model and our parameter 
estimation study to higher Post-Newtonian order. 

The fastest time scale for the system is the orbital time 
scale, which to leading order is given by Kepler's law: 



orb 



(1) 



where M = mi + m2 is the total mass of the system and 
a is the semi-major axis. Periastron advance enters at 
1 PN order, with a fractional advance per orbit of fc = 
3AI^ / , where fj, = mini2/M is the reduced mass of 
the system, L ~ aM(l — e^) is the orbital angular 
momentum and e is the eccentricity of the orbit. The 
ratio of the periastron precession timescale to the orbital 
timescale is given by 
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The orbital angular momentum L precesses with an- 
gular velocity Wpicc = Soff / |35[ , where 



Scff' = 21 
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The effective spin vector has magnitude S^h ~ since 
the magnitude of the individual spins is given by Si = 
XiiTT-i where the dimensionless spin parameter < Xi ^ 
1. The ratio of the precession time scale and the orbital 
timescale is given by 
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(4) 



We see that the precession of the orbital plane enters 
at 1.5 PN order relative to the orbital time scale. The 
precession of the orbital angular momentum due to spin- 
orbit coupling results in a modulation of the amplitude of 
the gravitational waveform at the solar system baryccn- 
ter. 
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FIG. 2: The time dependent gravitational waveform h+{t) at 
the solar system barycenter for a binary black hole system 
with eo = 0.3 (Source 1 in Table |ll| during the final year 
before merger. 



The rate at which the binary black hole system loses 
energy and angular momentum due to gravitational wave 
emission defines the decay time scale 
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and the ratio 
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The loss of energy and angular momentum results in the 
decay of the semi-major axis and the radial eccentricity 
of the system. Orbit avera ged expression for these decay 
rates can be found in Rcf. ^5|- The effects of spin-orbit 
induced precession of the orbital plane and the overall 
sweep up in frequency and amplitude as the system in- 
spirals over the course of a year are apparent in Figure [2] 
for a system with eo = 0.3 and other source parameters 
given in Table HIl for Source 1. 

At times well before the merger of the system M/a is 
small and we find that 
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We take advantage of this separation of the relevant 
time scales to make our waveform calculations more effi- 
cient. Since the decay and precession timescales are much 
longer than the orbital timescale we can start with a so- 
lution to the orbital equations of motion that neglects 
dissipation. This allows us to use Nyquist sampling of 
just a couple of samples per orbit for the quantities that 
depend on the dissipation and precession equations and 
saves considerable computational cost. Considering only 
times early in the evolution of the system we have a clean 
separation of time scales and the different processes can 
be treated differently in the calculations. This simplifies 
the problem by allowing an adiabatic treatment where 
the dissipation is assumed to be small over the course 



of an orbit and the eccentricity and semi-major axis are 
treated as constant while calculating an individual orbit. 

The approximations involved in exploiting the sepa- 
ration of timescales introduce small errors in the wave- 
forms. The largest of these comes from using the orbit- 
averaged spin precession equations, which neglects small 
periodic changes in the spin orientations that occur on 
the orbital timescale, but these changes are effectively 
2.5 FN order contributions to the higher harmonics. The 
averaging also introduces small errors in the long term 
secular evolution that scale as the ratio of the averag- 
ing time scale (in our case the orbital timescale) and the 
timescale of the terms being averaged (such as the spin 
precession timescale). These errors are multiplicative, 
and so represent higher FN order terms that can be dis- 
carded at 1.5 FN order. 

As the system approaches merger the various time 
scales become comparable and our waveform model 
breaks down. We adopt the termination condition 
27rM/orb = 0.01, which corresponds to an expansion pa- 
rameter M/a K, 0.05. Denoting the time when this con- 
dition is met as tg, and the orbital frequency at this time 
as fs, we taper the waveform smoothly to zero by multi- 
plying the waveform with a half-Hann filter: 
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with tH^ts- 3/(2/,). 

Our termination condition is conservative in terms of 
the signal to noise ratio (SNR) LISA will be able to ex- 
tract from this type of source. Most of the SNR comes 
from times near merger, so the extension of the validity 
of the waveform closer to merger results in a big increase 
in SNR. Our study here is thus a pessimistic estimate 
of how well LISA will be able to determine the various 
source parameters. In the case of the radial eccentricity 
parameter however, the circularization of the waveform 
toward merger means that most of the eccentricity infor- 
mation is encoded at times well before merger. While 
increased SNR would improve the determination of the 
other source parameters, we find that the eccentricity is 
not highly correlated with the other parameters (see Fig- 
ure [T3]). Our choice for when to truncate the waveform 
thus should not have much of an effect on our study of 
how well LISA will be able to determine the eccentricity 
of black hole binary systems. 

We have tested our waveform generator in various lim- 
its against other codes. In the circular limit, and with 
dissipation turned off, we found precise agreement with 
the 1.5 FN limit of the spinning black hole code developed 
by Cornish, Hughes, Lang and Nissanke that is described 
in Ref. [s^l ■ We do not expect, and nor do we find, precise 
agreement when dissipation is included. This is because 
our eccentric waveform generator evolves the semi-major 
axis, while the circular orbit code evolves the orbital fre- 
quency, which leads to numerical differences at 2 FN or- 
der. In the 0-FN limit (hence no spin effects and no 
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higher harmonics) wc find precise agreement with the 
Peters and Mattliews waveforms Q- 



III. LISA RESPONSE 

We simulate the LISA response to a gravitational wave 
signal plus instrument noise and confusion noise due to 
galactic binary sources of gravitational waves. We adopt 
the standard ecliptic coordinate system with origin at the 
barycenter. The individual data streams from the six 
LISA phase meters can be combined to cancel out the 
laser phase noise and form Time Delay Interferometry 
(TDI) variables The Michelson style TDI variables 
{X, Y, Z} can be used to construct three noise orthogonal 
data streams that are similar to the {A, E, T} variables 
described in Ref. [soj . 

Processing the gravitational waveform through the 
LISA response function imparts additional amplitude 
and frequency modulations on the one year timescale of 
the LISA orbits. These effects can be seen in Figure |31 
which shows the A-channcl response to the Barycenter 
signal shown previously in Figure [5] 
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FIG. 3: The A{t) channel response to the spinning binary 
black hole system shown in Figure [2] The amplitude modula- 
tion is due to a combination of the antenna pattern sweeping 
around as LISA orbits the Sun, and the spin induced pre- 
cession of the orbital plane. The overall gravitational wave 
amplitude grows as the system spirals in and nears merger. 

The one-sided noise spectral density of the detector in 
the A and E channels is given by [4ll |: 



5i„st(/) = -[{2 + cosu)Sm 
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where u = 27r//L, SfKf) is the position noise and S'^{f) 
is acceleration noise, and we have taken the limit of sym- 
metric noise in the detector. The confusion noise due 
to gravitational wave sources in the galaxy has been es- 
timated by direct simulation of the LISA response [i^ 



to a synthetic ga laxy [43 | , followed by the removal of re- 
solvable systems [4J] . An approximate fit to the resulting 
confusion noise estimate is given by 



5'conf(/) = < 



' 10~''4.8 y-2-4 
-[^Q-47.15 J-3-1 

10-74.7J-13 
2Q-59.15 j-7 



10-* < / < 4.5 X 10-4 
4.5 X 10-'* < / < 1.1 X 10-3 
1.1 X 10-3 < / < 1.7 X 10-3 
1.7 X 10-3 < / < 2.5 X 10-3 
2.5 X 10-3 < / < 4 X 10-3 
(10) 

The total noise is then taken to be the sum of the the 
two contributions: 5'„(/) = Sinstif) + Sconiif)- 

IV. PARAMETER ESTIMATION 

The goal of a parameter estimation study is to find 
how accurately we can determine the values of the source 
parameters for a given signal. We take a Bayesian ap- 
proach and use the Markov Chain Monte Carlo (MCMC) 
technique to compute the posterior distribution function 
p{x\s) describing the model parameters x that we infer 
from data s. Our results establish how well the various 
source parameters will be able to be determined by the 
LISA mission, including the orbital eccentricity. 

The use of MCMC techniques in gravitational wave 
data analysis is now a familiar technique for parameter 
estimation in gravitational wave astronomy |45l450| . The 
result of a well constructed MCMC is a set of samples 
from the posterior distribution. The number of samples 
from a particular region of parameter space is propor- 
tional to the posterior weight contained in that region. 
The uncertainty in each parameter is given by quantiles 
of the marginalized posterior distribution (e.g. the half- 
width of the 90*^ percentile equates to a 2-t7 error if the 
distributions are Gaussian). 

By Baycs theorem, the posterior distribution is given 
by the product of the prior distribution p(x) and 
likelihood p(s|x), normalized by the evidence p{s) ~ 
/p(x)p(x|s) dx. For Gaussian noise, the likelihood of the 
data s having being generated by a gravitational wave 
signal /i(x) is given by 



p(s|x) = Ce 



-i(s-/l(x)|s-?l(x)) 



(11) 



where C is a normalization constant that does not depend 
on the signal or the template. Here we have used the 
noise weighted inner product 



{a\b) 
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(12) 



where £'„(/) is the one-sided noise spectral density. The 
prior probability density for x reflects our knowledge of 
the source parameter, however ill-formed, before we an- 
alyze the data. For example, we assume a uniform prior 
for angular parameters such as the sky location and the 
initial orientation of the angular momentum vector, such 
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that the cosine of the co-latitude is uniformly distributed 
in the range [—1 : 1] and the azimuth is uniformly dis- 
tributed in the range [0 : 27r]. 

The primary mode of the posterior distribution yields 
the best fit values for the source parameters, according 
to the current data and our prior knowledge. In many 
instances the posterior distribution is multi-modal, and 
the quantiles used to estimate the parameter uncertain- 
ties may cover disjoint regions in parameter space. Even 
when the bulk of the posterior weight lies in a single, con- 
tiguous region, the posterior distribution may not be well 
approximated by a Gaussian distribution. Nonetheless, 
a Gaussian approximation to the posterior distribution 
often provides a reasonable estimate of the parameter es- 
timation errors, which can be efficiently computed using 
the Fisher information matrix Tij, which measures the 
expectation value of the curvature of the posterior distri- 
bution about the mode: 

= -{didj lnp(x|s))|,„ode • (13) 

We employ a parallel tempered [s^l Metropolis- 
Hastings [Hi [111 MCMC routine to explore the PDF. 
The Markov chain starts at parameter values x and tran- 
sitions to y with probability 



V. PARAMETER ESTIMATION WITH LISA 

We can use our time dependent gravitational wave- 
forms [11] and established MCMC techniques to study 
how well LISA will be able to measure the full set of 
seventeen parameters necessary to describe a spinning 
binary black hole system in an eccentric orbit. We are 
especially interested in determining when eccentric orbits 
can be distinguished from circular orbits. 

We chose to study signals in their final year prior to 
merger, with an observation time that extends just be- 
yond the merger. In order to choose initial parameter 
values for a system that will merge in one year, we need 
to calculate an initial semi-major axis based on the life- 
time estimate for a system with some given initial radial 
eccentricity. We use the leading order, OPN expression 
for ao i 

ff^M\Y ( 157 2 5799977 4 

Oo = 4 ic 1 H en H en 

" V 5 7 V 172 ° 7336832 ° 

1888175763 f\ 

H en , (15) 

2523870208 °) ^ ' 

as an initial guess, and apply a bisection routine to the 
full numerical orbital evolution to find the value of ao 
that yields a merger time of one year. 
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Here q(x.\y) is the proposal distribution, which is the 
function that generates proposals for moves from x to 
y. The performance of an MCMC algorithm is quite 
sensitive to the choice of proposal distribution, and care 
must be taken to ensure that the chains do not get stuck 
on local maxima of the PDF. We employed several tech- 
niques to ensure rapid exploration of the full parameter 
space: local coordinate transformations to uncouple the 
parameters; moves that exploit symmetries of the likeli- 
hood surface to encourage jumps between local maxima; 
and parallel tempering to encourage wide exploration of 
the posterior (50, . .51, .56| . The number of iterations spent 
at each parameter value is proportional to the posterior 
density, and histograms of the parameters visited by the 
chain provide an estimate of the posterior distribution. 

We have found that drawing from a variety of proposal 
distributions provides a set of jump proposals that tend 
to produce an MCMC that efficiently maps out the de- 
sired PDF and provides accurate parameter uncertainties 
even for very large search spaces. Our parameter estima- 
tion study thus uses several proposal distributions, in- 
cluding parallel tempering and Fisher matrix proposals. 
In this high dimensional parameter space we found it ad- 
vantageous to use the Fisher matrix approximation to the 
posterior to propose jumps along single eigen-directions 
as part of the mixture of jump proposals. 



TABLE L Parameter ranges for our study of spinning black 
hole binary systems in eccentric orbits. 



Parameter 


Minimum 


Maximum 


mi 




WMq 


m2 




IQ'Mq 


Dl 


1 Gpc 


100 Gpc 


eo 





1 


ao 


20 M 


1000 M 


Xi 





1 


X2 





1 


cos 9 


-1 


1 


cos 9l 


-1 


1 


cos Osi 


-1 


1 


cos 


-1 


1 







2tv 







2tv 







2tv 


4>S2 





2tv 


no 





2tv 







2tv 



We use parameter ranges consistent with typical LISA 
sources, given in Table HI The masses are given in terms 
of the mass of the Sun, Mq = 1.9891 x lO^o kg, and 
the luminosity distance Dj^ is given in units of Giga- 
parsecs. The dimensionless spin parameters xi ^-nd X2 
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combine with the black hole masses to give the magni- 
tudes of the spins, Si = Xi''^i- There are initial ori- 
entation parameters for the orbital angular momentum 
vector L — > {cos 9 l , 4> l) , as well as the spin vectors 
Si — > (cos 9 Si , 4'Si)- The final two parameters, no and (jjQ, 
are initial phases related to the mean motion and orbital 
phase. In the circular limit these parameters arc degener- 
ate, but for eccentric orbits we have to specify the initial 
periastron position. In the present study we have ignored 
the possibility that gas dynamics may partially aligned 
the spins with the orbital angular momentum |57| , which 
would restrict the prior ranges and reduce the degree of 
orbital precession. The impact of this partial alignment 
on parameter estimation has been considered for circu- 
lar orbits [EH , and it would be interesting to extend this 
study to include eccentricity. 

Here we study several representative cases to estab- 
lish the parameter recovery errors and to study corre- 
lations between the parameters. The high dimension of 
the parameter space makes it difficult to perform a com- 
prehensive study - if we were to choose just two values 
of each parameter we would need to perform 2^^ 10^ 
parameter estimation studies. Each MCMC run involves 
^ 100, 000 iterations with ^ 8 parallel chains, and takes 
about a week to run on a single 2.66 GHz Intel processor, 
so we are limited in the number of examples we can con- 
sider. We perform MCMC parameter estimation studies 
of several representative examples, varying the mass ra- 
tio, sky location, distance, eccentricity, and dimensionless 
spin parameters. We only looked at a few initial spin and 
orbital orientations since we do not expect these to have 
a significant effect on the results - unless the initial ori- 
entations are very special the system will explore a wide 
range of orientations during the orbital evolution. To be 
able to explore the parameter space more widely would 
take a faster code. One possibility is to use the Fisher in- 
formation matrix approximation to the posterior, which 
is many orders of magnitude faster than a full MCMC 
study. In preparation for such a study we compare the 
Fisher matrix approximation to the MCMC derived pos- 
terior distributions and find that the approximation is 
fairly reliable so long as the initial eccentricity exceeds 
eo - 0.01. 

Our first study focuses on determining when the ec- 
centricity one year before merger is distinguishable from 
zero. We studied two systems that only differed in sky 
location (and hence in SNR), and considered initial ec- 
centricities in the range eo S [0.001,0.2], see Table HIl for 
a list of the other source parameters. Marginalized poste- 
rior distributions for cq are shown in Figure H] for Source 
1. and Figure [5] for Source 2. We see that the error in 
the measured value of eo gets smaller as eo gets larger, 
but the dependence on eo is fairly weak, and never larger 
than Aeo ~ 0.001. Our criteria for deciding if the eccen- 
tricity is distinguishable from zero is to see if their is any 
weight in the posterior distribution at eo = (this test is 
motivated by the Savage-Dicke density ratio estimate for 
the model evidence 43 1) . For Source 1 we see that the ex- 



TABLE II: Injected parameter values for two sets of sources 
studied with a range of vahies for gq. The results of the pa- 
rameter estimation study for Source 1 are given in Figure |4] 
and the results for Source 2 are given in Figure [S] 



Parameter 


Source 1 


Source 2 


mi 


2 X 10^ Mq 


2 X 10^ Mq 


m2 


1 X 10^ Mq 


1 X 10'' M0 


Dl 


6.36167 Gpc 


6.36167 Gpc 


Xi 


0.5 


0.5 


X2 


0.8 


0.8 


COS 6 


0.2 


0.4 


cos 6l 


-0.5 


-0.5 


cos 


-0.8 


-0.8 


cos 6 $2 


0.6 


0.6 


<!> 


1.2 


2.0 




2.6 


2.6 




0.4 


0.4 




1.7 


1.7 


no 


0.2 


0.2 


00 


1.65 


1.65 



amples with eo > 0.005 are clearly distinguishable from 
circular, while the eo = 0.002 case is on the margin of de- 
tectability. For Source 2, which has higher SNR due to a 
more favorable sky location, the eo = 0.002 case is clearly 
distinguishable from circular. The Fisher matrix approx- 
imation to the posterior distribution is computed at the 
maximum a posteriori probability (MAP) value of the 
parameters, and is found to work well for eccentricities 
eo > 0.05, but breaks down for small eccentricities. The 
MCMC derived posterior distributions are much fiatter 
than a Gaussian distribution, and we attribute the fail- 
ure of our Fisher matrix estimates to only including the 
leading order, quadratic curvature terms in the Fisher 
matrix calculation. 

This study suggests that LISA observations of eccen- 
tric black hole binary systems will be able to measure 
the eccentricity of the system and distinguish eccentric 
systems from circular systems to parts in a thousand. A 
similar study was performed by Porter and Sesana [l^l for 
non-spinning eccentric binary black hole systems. Their 
results suggest that LISA will be able to measure the 
eccentricity to parts in 10~^ for such sources. 

We find that the other source parameters arc also mea- 
sured quite well, as illustrated in Figure [5] and Figure [T] 
Marginalized posterior distributions are shown for the 
chirp mass Mc = (mim2)'^/^/(mi -I- 7712)^^^ and reduced 
mass ^ = toiTO2/(toi + ^2), the distance to the source, 
the initial radial eccentricity, and the two sky location 
parameters. We compare the Fisher matrix approxima- 
tion to the marginalized posterior distributions computed 
from the MCMC runs and find excellent agreement for all 
the parameters (except for the eccentricity in Figure [B]) . 
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0.009 0.01 0.011 0.012 



0.002 0.004 0.006 0.008 0.001 0.002 0.003 

FIG. 4: Marginalized posterior distribution for the initial ra- 
dial eccentricity for sources with the same parameter values 
except for different initial eccentricities and semi-major axes. 
The boxed histograms are derived from the Markov chains, 
while the solid lines are the Fisher matrix predictions com- 
puted at the MAP values of the parameters (which are pushed 
off the true values by the instrument noise). Top left eo = 0.1, 
ao = 69.96M, SNR= 231; top right eo = 0.01, ao = 68.33M, 
SNR= 237; bottom left eo = 0.005, ao = 69.3M, SNR= 237; 
bottom right eo = 0.002, ao = 69.3M, SNR= 237. The other 
parameter values correspond to Source 1 in Table HIl 



0. 002 0.1004 




FIG. 5: The marginalized posterior distribution for the initial 
radial eccentricity for sources with the same parameter values 
except for different initial eccentricities and semi-major ax;es. 
On the left eo = 0.1, ao = 68.94M, SNR= 557; on the right 
eo = 0.002, ao = 68.32Af, SNR= 559. For this source, eo = 
0.002 is distinguishable from the circular case. A vertical line 
marks the injected eo value in each plot. The other parameter 
values are given in Table IHl Source 2. 



As a start to exploring the large parameter space of 
eccentric binary black hole systems we consider a few 
representative examples below. We find that in general 
physical parameters are well constrained by LISA obser- 
vations and that the Fisher matrix makes fair estimates of 
the parameter errors. In addition to the examples below 
with varied spin and distance parameters we also studied 
systems with a range of mass ratios mi/m2 £ [1,5] and 
total masses M S [IO^Mq, IO^Mq] with similar results. 
We also compared the parameter estimation errors with 
those obtained when the eccentricity is held fixed at zero, 
and saw only small (less than 50%) changes in the error 
estimates. 

The dimensionlcss spin parameters xi ^-nd X2 a-re var- 
ied in Table Hill to study the cases of large and small spin 




14.0116 14.0116 14.0117 14.0117 





0.0001469890.00139639 0.00264579 0.00389519 





0.192008 0.200016 0.208024 0.216031 



FIG. 6: The marginalized posterior distribution for several 
source parameters, including the initial radial eccentricity 
with injected value eo = 0.002, ao = 69.3M, and SNR= 237. 
The other injected source parameters are those of Source 1 
in Table HTl The solid lines are the Fisher matrix predictions 
computed at the MAP values of the parameters. 





14.0116 14.0116 14.0117 14,0111 
ln(DL/Gpc) 



13,4098 13.41 13.4101 13.4103 





1.84227 1.85383 

00B(8) 



0.0995662 0.0998871 0.100208 0.100529 





0.381348 0.392921 0.404493 0.416066 



FIG. 7: The marginalized posterior distribution for several 
source parameters, including the initial radial eccentricity 
with injected value eo = 0.1, ao = 68.94M, and SNR = 558. 
The other injected source parameters are those of Source 2 
in Table [nl The solid lines are the Fisher matrix predictions 
computed at the MAP values of the parameters. 



parameters. The magnitude of the spin vectors is related 
to the mass of the black hole Si = Ximf. The poste- 
rior distribution for several of the source parameters for 
these cases are given in Figures [5] and [HI The Fisher ap- 
proximation is a reasonable prediction of the width of 
the posterior distribution for the spin of the more mas- 
sive body, but does a poor job for the less massive body. 
This discrepancy was seen in many other examples that 
we looked at. The cause of the discrepancy is presently 
not understood. Two lower SNR examples are shown in 
Figurc[Tn]and FigurelHl and we see good agreement with 
the Fisher matrix estimates. 

We do not expect the eccentricity to be highly corre- 
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TABLE III: Source 3 gives the injected parameter values for 
Figure [S] with large spin values. Source 4 gives the injected 
parameter values for Figure |9] with small spin values. 



Parameter 


Source 3 


Source 4 


nil 


2 X IO^Mq 


2 X 10^ Mq 


m2 


1 X IO^Mq 


1 X lO^'M© 


Dl 


6.36167 Gpc 


6.36167 Gpc 


eo 


0.1 


0.1 


ao 


68.96 M 


68.91 M 


Xi 


0.8 


0.1 


X2 


0.9 


0.11 


cos 6 


0.2 


0.2 


cos 9l 


-0.5 


-0.5 


cos 


-0.8 


-0.8 


cos 


0.6 


0.6 




1.2 


1.2 




2.6 


2.6 




0.4 


0.4 




1.7 


1.7 


no 


0.2 


0.2 




1.65 


1.65 


SNR 


250 


197 


In(M^) 




ln(,i) 







0.795638 0.802617 0.809595 0.816574 

oos(e) 



0.867984 0.895003 0.922022 0.949041 





0.192361 0.197012 0.20161 





14.0118 14.0119 




— 1^ — Tirh^^^ 



0.0755431 0.0879182 0.100293 0.1 



0.0993269 0.105033 0.110739 0.116444 





0.166536 0.187678 0.20882 0.229961 



FIG. 9: The marginalized posterior distribution for several 
source parameters for Source 4 in Table lXITl The solid lines are 
the Fisher matrix predictions computed at the MAP values 
of the parameters. 



TABLE IV: Source 5 gives the injected parameter values for 
Figure [1171 with redshift z ~ 1.5. Source 6 gives the injected 
parameter values for Figure [11] with redshift z ^ 2. 



FIG. 8: The marginalized posterior distribution for several 
source parameters for Source 3 in Table [Till The solid lines are 
the Fisher matrix predictions computed at the MAP values 
of the parameters. 



Parameter 


Source 5 


Source 6 


mi 


2 X 10*^ M© 


2 X 10" Mq 


m2 


1 X 10" M© 


1 X 10" Mq 


Dl 


11.008 Gpc 


15.733 Gpc 


eo 


0.1 


0.1 


ao 


68.94 M 


68.94 M 


XI 


0.5 


0.8 


X2 


0.5 


0.8 


cos 9 


0.2 


0.2 


cos Ol 


-0.5 


-0.5 


cos 9si 


-0.8 


-0.8 


cos 6^52 


0.6 


0.6 




1.2 


1.2 




2.6 


2.6 


<t>Si 


0.4 


0.4 


<i>s2 


1.7 


1.7 


no 


0.2 


0.2 


4>o 


1.65 


1.65 


SNR 


132 


92 



latcd with the other source parameters since the higher 
harmonics introduced in the waveform due to eccentric- 
ity cannot be simulated by changes in other parameters 
or their combinations. We indeed find that the initial 
eccentricity is not correlated with the other parameters. 
Compare the distribution of values for the two masses 
in Figure [T^ to the distribution of mass and eccentricity 
values in Figure [T51 The two masses are highly correlated 
since it is the total mass of the system and the ratio of 



the masses that appear in the waveform. The distribu- 
tion of eccentricity versus the other source parameters is 
similar to that seen in Figure 1131 



VI. CONCLUSION 

Our studies of the response of the LISA detector to 
the gravitational wave signal from spinning binary black 
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14.0114 14.0116 14.0119 14.0121 




2.37793 2.40825 2.43857 2.4 




0.189715 0.203912 0.218109 0.232306 





0.0979112 0. 



0.100041 0.101107 




FIG. 10: The marginalized posterior distribution for several 
parameters for a source with z ~ 1.5 eo = 0.1, and SNR = 132 
f Table HVl Source 5). The solid lines are the Fisher matrix 
predictions computed at the MAP values of the parameters. 




14^ 14.SCK 14^09 14.£1 14.^11 



14.aJ7 14.SQB 14.509 14.51 14.511 
h(ml) 



S a - 



aS4 ase osb ae o b2 o.64 o.es 




FIG. 12: Two dimensional posterior distribution scatter plots 
showing the correlation between pairs of parameters for a 
source with initial radial eccentricity eo = 0.002 and SNR 
= 237 (Table [ill Source 1). Upper left mi and m2, upper 
right mi and ao, bottom left cos^s^ 4^ bottom right 
no and 0o. 




0114 14.0117 14.012 




2.72976 2.76318 2.7966 
cos(O) 




0.162017 0.184954 0.207892 0.230829 




13.409 13.4098 13.4105 




0.0984195 0.0998153 0.101211 0.102607 




1.17681 1.19042 1.20403 1.21764 



FIG. 11: The marginalized posterior distribution for several 
parameters for a source with 2 ~ 2, eo = 0.1, and SNR = 92 
(Table HVl Source 6). The solid lines are the Fisher matrix 
predictions computed at the MAP values of the parameters. 



f aoaa 
aopi 



Q.0D3 

aocE 
aoni 




14Sa7 14.SaB 14.SDg 14.SI 14.511 



Q2 a2i 




o.oce 
aoQi 



6928 egjai eaaa eg.35 63.37 



0.47 0.4S CL« as 0-61 OSS 



FIG. 13: Two dimensional posterior distribution scatter plots 
for a source with initial radial eccentricity eo = 0.002 and 
SNR = 237 (Table im Source 1). Upper left mi and eo, upper 
right ao and eo, bottom left cos 9 and eo, bottom right xi and 
eo- 



hole systems in eccentric orbits show that the eccentric- 
ity should not be neglected for LISA data analysis and 
parameter estimation. We find that LISA can determine 
the eccentricity of the system one year before merger to 
parts in a thousand. This result depends only weakly on 
the initial value of the eccentricity, indicating that LISA 
will be able to distinguish between eccentric and circular 
orbits at the same level (eo ^ 10"'^). 

The construction of the gravitational waveforms for 
spinning binary black hole systems in eccentric orbits to 
1.5 PN order in [s^ establishes the framework for the ex- 
tension of this work to higher post-Newtonian order. Bi- 
nary black hole waveforms that include eccentricity will 
be necessary for several of the LISA science goals in- 
cluding constraining galaxy merger scenarios and testing 
general relativity in the strong field regime near super- 
massive black holes. 



This work builds the foundation for further studies 
including a comprehensive exploration of the parame- 
ter space including source sky location, distance, spins, 
masses, and mass ratios. For sources with initial radial 
eccentricity greater than cq ~ 0.02 the Fisher matrix is 
a good approximation to the posterior distribution. The 
very large parameter space could be studied quickly us- 
ing the Fisher approximation, although it is not as useful 
for the low initial eccentricity cases. 

The equations of motion and instantaneous gravita- 
tional waveforms have been calculated to the next post- 
Newtonian order and these pieces can be included in a 
future parameter estimation study. The 2PN effects in- 
clude the spin-spin coupling of the two black holes and 
thus corrections to the precession and evolution of the 
system. 

These waveforms can be also be used to study how well 
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the Advanced LIGO- Virgo network and proposed Ein- 
stein Telescope will be able to measure eccentricity and 
what level of bias could be expected from using circular 
templates for parameter estimation. 



to Alberto Scsana and Edward Porter for extensive and 
informative discussion of the mechanisms responsible for 
producing black hole binaries with eccentric orbits. 
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